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Abstract 



The disordered Bose Hubbard model is studied numericahy within the 
Bogohubov approximation. First, the spatiahy varying condensate wavefunc- 
tion in the presence of disorder is found by solving a nonlinear Schrodinger 
equation. Using the Bogoliubov approximation to find the excitations above 
this condensate, we calculate the condensate fraction, superfluid density, and 
density of states for a two-dimensional disordered system. These results are 
compared with experiments done with ^He adsorbed in porous media. 
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I. INTRODUCTION 



By definition, superfluids are robust to the introduction of weak microscopic disorder. 
A flowing superfiuid is cliaracterized by a macroscopic wavefunction whose phase varies 
across the sample; as long as this condensate wavefunction remains well defined, disorder 
cannot lead to the degradation of currents for topological reasons.0 With increasing disorder, 
however, the rigidity of the superfiuid towards phase variations is reduced. For sufficiently 
large disorder superfluidity is eventually destroyed even at zero temperature, resulting in a 
Bose insulator (the "Bose glass" 

Disordered Bose condensates can be realized experimentally by superfluidsB in random 
media, such as "^He films adsorbed on porous Vycor glass. The first few monolayers of 
adsorbed ^He are not superfiuid, even at low temperatures, and form an "inert" insulating 
layer of bosons localized by disorder. As the coverage is increased, a transition from this 
Bose insulator to a superfiuid phase is observed.! Crudely speaking, the first few monolayers 
are comprised of bosons occupying non-overlapping localized states, which screen the mi- 
croscopic disorder of the porous glass for subsequently added bosons. Added bosons feel a 
smoother potential that is the sum of the initial random potential plus a Hartree repulsion 
from the localized particles. When the disorder is sufficiently well screened, condensation 
into an extended state occurs.0 

Of course, this picture is an oversimplification: the ^He atoms in the "inert" layer are 
indistinguishable from those in the condensate, and the true many-body wavefunctions must 
be completely symmetric with respect to particle interchange. Exchange between the "inert" 
and "condensed" bosons can be important, especially near the insulator-superfluid transi- 
tion. The computational problem with this scheme is that, unlike the case of fermions, 
which by the exclusion principle must populate orthogonal states, bosons actually prefer to 
be in non-orthogonal states to optimize their effectively attractive exchange interactions. 
The need to symmetrize thwarts controlled general Hartree-Fock calculations. 

We present here numerical calculations of the properties of highly disordered Bose 
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condensates using the Bogoliubovllj approximation, which has been formulated for disor- 
dered systems by Lee and Gunn0 and considered in the weak disorder hmit by Meng and 
Huang.Ei'0 Ahhough the Bogohubov method is strictly valid only in the limit of weak re- 
pulsive interactions, we will consider strongly interacting systems as well (in the presence of 
arbitrary disorder) in an attempt to address the qualitative features of Bose systems in ran- 
dom media in a quasi- analytic fashion. Previous theoretical approaches include numerical 
simulations, 0^ scaling analysisjl renormalization group calculations,i'0'@ and perturba- 



23 



tive methods-t 

In the Bogoliubov approximationlll the disordered potential is screened by bosons occu- 
pying a delocalized condensate wavefunction which has larger amplitude where the random 
potential is deep. This non-uniform condensate is macroscopically occupied, and fluctua- 
tions into and out of it are considered due to residual interactions. These effects deplete the 
condensate non-uniformly, and lead to a spectrum of collective, phonon-like excitations. 

The Bogoliubov scenario resembles the heuristic "inert layer" picture discussed above, 
but is constructed in the reverse order. First the condensate is determined, and then the 
(possibly localized) non-condensate part of the many-body wavefunction is considered. This 
localized, uncondensed part of the ground state corresponds to the "inert layer" discussed 
above, and can be crudely thought of as the zero temperature "normal" fluid excited from 
the condensate by disorder rather than thermal fluctuations (see Sec. VI for a more precise 
definition). The advantage of the Bogoliubov approach is that exchange between the "con- 
densate" and the "normal" fluid is included naturally. Its disadvantage is that interactions 
within the normal fluid are essentially ignored. 

The outline of this paper is as follows: In Section II, we introduce the disordered Hubbard 
model for bosons and in Section III we solve this model in the Hartree approximation. In 
Section IV we review the Bogoliubov approximation for disordered bosons. Sections V 
and VI present calculations of the depletion of the condensate and the reduction of the 
superfluid density due to disorder, respectively. Section VII reports calculations of the 
excitation spectrum and specific heat of the disordered condensate. Finally, in Section VIII 



we summarize our results and discuss experiments. 



II. THE BOSON HUBBARD MODEL 

A simple model for disordered interacting bosons is the Hubbard model for lattice bosons 
in a random potential: 

■H = -tY: b\b, + y: v{i)b\h + y E ^1^1^^^- (1) 

where h\ {pi) creates (destroys) a boson at lattice site i. The sum on extends over all 
nearest-neighbor pairs of lattice sites, U is the strength of the repulsive on-site interaction, 
and t is a hopping matrix element. The random potential V{i) is uniformly distributed 
between —A and A. The total number of bosons is A/", and the number of lattice sites is V; 
the mean density is then n = M /V. 

As a model for the behavior of ^He adsorbed in Vycor on length scales less than the pore 
size (several hundred Angstroms), each site could represent a surface location of atomic 
dimension, connected to neighboring sites in a two-dimensional network. We will consider 
(HI) on two-dimensional square lattices of up to 306 sites, with periodic boundary conditions. 
To study disordered Bose condensates on longer length scales, the sites of model (|l|) could 
themselves be used to represent pores in Vycor, with a three-dimensional connectivity. Since 
300 sites is still quite a small three-dimensional lattice, we will only report calculations in 
two dimensions. 

III. THE HARTREE CONDENSATE 

A simple variational ground state for (|l]) is the Hartreeii state 

^'(ri,r2, ...,r7v) = 0o(n)0o('^2)---0o('^7v), (2) 

where all bosons are condensed into the same real, normalized single-particle wavefunction 
00 (0- The many-body state (H) is explicitly symmetric under particle exchange, as befits 
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a Bose state. For a translationally invariant system, the single-particle state 4>o{i) is inde- 
pendent of i, and is simply the zero-momentum state. In a disordered system, this will no 
longer be the case: (f)o{i) will adjust to be larger at the minima of the random potential and 
smaller at its maxima. 

The expectation value of the Hamiltonian (|I|) in the variational state is 

+ ^E'^o(^)■ (3) 

i 

To minimize (|^) with respect to the (normalized) single particle state 0o one must solve the 
discrete nonlinear Schrodinger equation 

-t E Mj) + vit)M^ = if^o + ex)M^), (4) 

j=nn(j) 

where the sum over j = nn{i) extends over the nearest neighbors j of site i. The effective 
single-particle potential V{i) is given by 

V{z) = V{z) + UN\M^)\', (5) 

where 4>o{i) is the single particle ground state of (^). 

For convenience, /io in eq. @ is chosen so that eo = 0, i.e., so that the Hartree excitation 
energies ex are measured with respect to the energy required to add a particle to the con- 
densate. The condensate 0o(^) and the V — 1 excited states denoted by (t>\{i) together form 
an orthonormal basis for single particle states. For convenience, sums over A will always 
implicitly exclude the condensate. 

We solve (^) and iteratively, as follows. Beginning with a trial condensate (f)o{i) (ei- 
ther the zero-momentum state or the exact state for t = 0), we compute the corresponding 
screened potential V{i). The resulting single-particle Schrodinger equation is solved numer- 
ically to obtain a new set of single particle eigenstates. An improved trial condensate is 
then created by mixing the initial guess with the lowest energy eigenstate of the screened 
potential. This procedure is repeated until (^ and (H) are simultaneously satisfied. Simple 
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linear interpolation to obtain a new trial condensate converges very slowly, if at all. More 
rapid, consistent convergence was obtained with the Broyden mixing method commonly 
used in electronic structure calculations.@ Achieving convergence is the biggest obstacle in 
our calculation, particularly for large disorder, and limits the system sizes we can consider. 

The condensate wavefunction (j)o{i) accumulates at the minima of the applied random 
potential, so that the screened potential V{i) is smoother than V{i), with shallower minima. 
These minima of the screened potential are more uniform than those of the original potential, 
with approximately the same depth. Roughly speaking, variations in the condensate conspire 
to create a screened potential which resembles the initial random potential, but with its 
deepest minima lopped off, as shown in Fig. |l]. As the density n (or the interaction strength 
U) is increased, the minima of V{i) become shallower and shallower, since the screening is 
then more efficient. No long-range correlations are introduced in the screening process, as 
shown in Fig. |^. 

To address the nature of the condensate and Hartree excited states, we calculate the 
participation ratio 

of state A, which measures the number of sites at which (f)x{i) is appreciable. 

The condensate wavefunction (f)o{i) is always extended,! and "participates" in a finite 
fraction of the lattice sites. The extended nature of the self-consistent ground state of 
is demanded by the following argument: Assume that the state (f)o{i) were localized. Then 
for a non-vanishing density of bosons, macroscopically occupying this single-particle state 
as in (0) would confine a macroscopic number of interacting particles to a finite volume 
(the localization volume of </>o(0)- The interaction energy of the resulting many-body state 
would then vary as the square of the total particle number. In the thermodynamic limit, 
however, the total energy should be extensive. Thus the assumption of a locahzed 0o(O 
yields a contradiction, and the condensate must be extended. 

This argument does not preclude a condensate wavefunction (f)o{i) which is "lumpy" - 
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i.e., one which is a (nodeless) superposition of well-separated, localized states. Strictly 
speaking, such a lumpy state is extended and the resulting Hartree state is still a Bose 
condensate. We will see below that these lumpy condensates (found for strong disorder 
and weak repulsion) are particularly susceptible to depletion from scattering out of the 
condensate, and have substantially reduced superfluid densities. An alternative, and perhaps 
better, variational state could be constructed by instead placing a few particles in each of 
a large number of localized states. As these localized states overlap (they need not be 
orthogonal), however, it becomes difficult to calculate the energy and other properties of 
the properly symmetrized state. Unfortunately, the numerical tricks which enable efficient 
calculation with determinants fail for permanents. 

It is well known that all eigenstates of a generic random potential in one and two dimen- 
sions are localized. How then can the condensate 0o always be extended? The loophole that 
permits this is that the screened potential V{i) is not generic, but has been tailored to the 
problem at hand specifically to produce an extended ground state. The condensate is not 
a "typical" state, but one whose peaks and valleys have been fed back into the disordered 
potential itself via (^). The extended nature of the condensate does not violate any accepted 
lore of localization. 

An analysis of the participation ratio for (two-dimensional) systems ranging from V = 72 
to V = 306 suggests that for small disorder, all states are extended (i.e., have a localization 
length larger than our largest system). Fig. ^ shows that the participation ratio scales 
with the size of the system. For sufficiently strong disorder, we find that the participation 
ratios of the Hartree excited states become independent of system size, indicating that they 
have all become localized. Only the condensate remains extended. It is interesting that 
the use of the self-consistent potential V{i) converts the state of lowest energy (which in a 
typical single-particle localization problem would be the most localized state) to the unique 
extended state. 
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IV. THE BOGOLIUBOV APPROXIMATION 



Given the self-consistent Hartree condensate, (poli), we can proceed with the Bogohubov 
approximation. Following Lee and Gunn,0 we expand the boson field operator 6, in the 
complete set of operators 60 and {bx}'- 

&. = 0o(^)feo + E<^AW&A. (7) 
A 

Although interactions and the disordered potential will both deplete the condensate, in the 
Bogoliubov approximation this depletion is assumed to be small enough that the single- 
particle state 4>o{i) is still occupied by a macroscopic number A/q bosons. To order I/Aq we 
can then replace the creation and annihilation operators 6q and 60 for this state by vCVq. 
The total number of bosons in the system is the sum of those in the condensate and those 
not in the condensate: 

Af = Afo + J2bih. (8) 

A 

Expanding to first order in the depletion of the condensate, we then find 

6o = v^=VAr-^E b{b, + ... (9) 



Inserting (|^) and (||) into the disordered Hubbard model (|1]), and retaining all terms 
second order in b\ and 6a, yieldS 

= AT/Xo - ^ E \M^) r + E ^A&I&A (10) 
^ i XX' 

+ ^ E Sxx'ibibx' + bibx + bibl + bxby). 

^ XX' 

The first line of (|1^) specifies the single-particle and self-interaction energies of the conden- 
sate, and the energy t\ for adding a particle in the excited state A. The second line involves 
the inner product 

^AA'=E<^aW|0o(^)P0A'(^) (11) 



of states A and A' weighted by the condensate density, which gives the amphtude for (a) 
single-particle scattering by the condensate and (b) pair scattering into and out of the 
condensate. Since the condensate is non-uniform, these scattering processes will generally 
not conserve momentum. 

To arrive at (^) terms cubic and higher order in field operators bx and b\ have been 
discarded. This is equivalent to the random phase approximation, and includes interactions 
between the non-condensate bosons and the condensate while neglecting interactions among 
the uncondensed bosons. These approximations are controlled in the dilute, weakly inter- 
acting limit in which the condensate fraction Mq/M is close to unity (see Sec. V). Here we 
will push the Bogoliubov approximation to its limits, and hope that the qualitative results 
are representative of disordered Bose condensates. 

The quadratic Hamiltonian ([ToD can be diagonalized by canonical transformation to a 
set of quasiparticle creation and annihilation operators 7^ and 7 such that 



[^B,7i] = ^M7i, 



(12) 



where uo^ is the quasiparticle excitation energy. This transformation is accomplished by 
taking linear combinations of creation and annihilation operators: 



7, 



(13) 



Like the index A labeling the Hartree states, the index /x labeling quasiparticle states runs 
from 1 to V — 1. 

To satisfy (|T2|), the coefficients u^x, v^x must obey the generalized eigenvalue equation 



^XX' -DXX' U^x' 

J \ J 



-Bxx' Axx' 
( 



UJ,, 



Sxx' 
-6xx' 



U^j^X' ^ 



(14) 



where 



^AA' = Ca^AA' + UMSxX', 



(15) 



(Summation over the repeated index A' is implied.) Note that if (") is a solution with 
excitation energy u (corresponding to 7^), then (^) is a solution with —u (corresponding to 
7). The orthonormality conditions 

/ 

dwi (I 1I...J\J 

(16) 



u^x v^x 



\ 



J 



Sxy 
-6xx' 

are automatically satisfied by normalized solutions of (p^, and guarantee that the quasipar 
tide operators 7^^ obey Bose commutation relations: [7^, 7^^/] = 6^^' and [7^, 7^/] = 0. 
The ground state energy Eq in the Bogoliubov approximation is 



Eg = Ni2c 



■E 



(17) 



+ E [(^AA' + UNSxx')Vf^xv^,x' - UNSxx'U^xv^^x' 

XX'/j. 

The last line gives the zero-point contribution of the quasiparticle modes. 



V. THE CONDENSATE FRACTION 

The ground state wavefunction \G) in the Bogoliubov approximation is the state anni- 
hilated by all of the quasiparticle destruction operators 7^: 



\G) = {bl)^mxx'exp[-Mxx'bib{,]\vac), 
where |vac) is the state with no bosons and Mxy is defined implicitly by 

J2ufi\Mxx' = v^x'- 

X 



(18) 



(19) 



The number of particles in the condensate, A/q, is determined by calculating the mean 
number of bosons not in the condensate {J2x ^a^a) and subtracting it from the total particle 
number (see eq. (||)). 
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In a translationally invariant system, the condensate fraction measures the occupation of 
the zero-momentum state. In a disordered system, the proper definition of the condensate 
fraction is the largest eigenvalue of the one particle density matrix.ii The condensate density 
is then given by the square of the off-diagonal long-range order parameter: 

hm {G\blb,\G) = {G\bl\G){G\b,\G) 

= AfoM^)Mj)- (20) 

Fig. ^ shows the condensate fraction, Ao/A/", as a function of disorder for several values of 
the interaction strength. The calculation is done by averaging over 7 realizations of disorder 
on Lj. X Ly lattices where L ranges from 8 to 18. We then extrapolate to the thermodynamic 
limit. Even in the absence of disorder, particles are scattered out of the condensate as a 
result of their mutual interactions, and Mq/M is less than unity. For weak disorder, the 
number of bosons in the condensate stays roughly fixed, while the condensate wavefunction 
itself is distorted to accommodate the random potential. 

This insensitivity of the condensate fraction to weak disorder is a crude criterion for 
superfluidity, although the proper quantity to consider is the superfiuid density (see below). 
Fig. ^ shows that as the interaction strength Un/t increases, the system becomes more 
robust to the addition of disorder, so that larger values of A are required to further deplete 
the condensate beyond the effect of interactions alone. 

At large values of disorder, the condensate fraction drops to zero. As the condensate 
fraction becomes small, the approximation of truncating the Bogoliubov Hamiltonian (p!OD 
at quadratic order becomes worse and worse, and our calculations cannot be considered 
quantitative. Nevertheless, our calculation suggests that for sufficiently high disorder the 
condensate is destroyed, and a "Bose glass" is reached. Although the logic leading to it 
breaks down for A/q = 0, the Bogoliubov ground state (|T8|) with vanishing condensate 
density is a potentially useful variational state for the Bose glass. 

In the Hartree calculation of Sec. II, the applied random potential V{i) is screened by 
A/'f/|0o(OP, which is equivalent to assuming that all of the bosons are in the condensate. 
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As the condensate is depleted, does this estimate of the screened potential continue to hold? 
To check this, we compare the ground state expectation value of the boson density at site 
i, {G\blbi\G), with the density obtained in the Hartree approximation, A/'|0o(OP- Fig. 
^ shows, the density in the Bogoliubov approximation faithfully tracks the density in the 
Hartree approximation. Even when the condensate is significantly depleted, the particles 
scattered from it remain in their original vicinity, and continue to screen the initial random 
potential as if they had remained in the condensate. 



The superfluid density of a Bose condensate distinguishes between the low-frequency, 
long-wavelength transverse and longitudinal responses of the system. (This quantity should 
not be confused with the condensate fraction discussed above, which is a ground state expec- 
tation value that measures the degree of off-diagonal long-range order.) A longitudinal probe 
corresponds to boosting the system, and the entire fluid responds. A low-frequency trans- 
verse probe corresponds to a slow rotation of the system, which only couples to the normal 
fluid, leaving the superfluid untouched. The superfluid density is defined simply as the dif- 
ference between the longitudinal and transverse response. In principle, a Bose system can be 
superfluid without possessing true off-diagonal long-range order, the canonical example be- 
ing the two-dimensional Bose liquid at non-zero temperature below the Kosterlitz-Thouless 
transition, which has only algebraic correlations. 

The zero-temperature, zero-frequency, current- current response tensor Xjj(q, = 0) is 
given by the Kubo formula@ 



The sum extends over all intermediate excited states m, and the lattice current operator 
J(q) is defined by 



VI. SUPERFLUID DENSITY 




(21) 




k 



(22) 
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where 6k = J^i^^^^'^i- In principle the direct evaluation of Xij is straightforward given a 
complete knowledge of the excited states |m). 

In the continuum, the longitudinal response in the long wavelength limit is required by 
the the /-sum rule to satisfy 

limxx.(gx) = -2tAA, (23) 

i.e., the entire fluid participates in longitudinal flow. On a lattice, the /-sum rule is 
modifiedii so that 

lim x.x(gx) = -2tKs = -2t[Af - J2 I7(^k)], (24) 

where ek is the tight-binding dispersion given below in (|30D. (This sum rule holds in the 
presence of arbitrary disorder.) Eq. ( p4D implicitly defines A/'eg. Even though M^s 7^ A/", the 
longitudinal response of a lattice system still corresponds to the entire fluid. 

The transverse response of a Bose liquid is only due to the "normal fluid," since the 
superfluid component of the system can only participate in irrotational [i.e., longitudinal) 
flow. Thus we can define the number of bosons in the normal fluid, A/'„, by 

limxxx(gy) = -2tUn. (25) 

q^O 

This definition also holds in the continuum, with t replaced by fi^ /2m. 

In a non-superfluid system, the long-wavelength longitudinal and transverse responses 
are identical. Superfluidity occurs when the two responses become different. The superfluid 
(number) density Ug is then defined by the difference between the longitudinal and transverse 
response, per unit volume: 

Ms Ks - Mn 



(26) 



V V 

The longitudinal and transverse response functions Xxx(q'x) and Xxx^Qj) can be easily 
calculated numerically in the Bogoliubov approximation. The excited states |m) entering 



(pi|) are then all one and two quasiparticle states. Because the Bogoliubov approximation 
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does not conserve particle number, the lattice /-sum rule ( p4D is not satisfied. Explicit eval- 
uation of ( pTD shows that limg^oo XxxiQ'^) tends to the number of bosons in the condensate, 
Ao, rather than the total particle number.^ (Note that Ao is not the same as Af^s in (p^.) 

The normal fluid density, obtained by explicit calculation in the Bogoliubov approxi- 
mation of the transverse response function Xxxilf), appears to be more trustworthy. A/"n 
vanishes in the translationally invariant case {V{i) = 0), as expected. We therefore fol- 
low ref. ^ and adopt ( p6|) as our operational definition of the superfiuid density in the 
Bogoliubov approximation. 

Fig. D shows the transverse response function for a 12 x 13 lattice averaged over 6 
realizations of various weak disorder. Such small systems were used because the evaluation 
of the response tensor Xy(q, c<j = 0) in the presence of disorder requires four nested sums 
over quasiparticle states for each q and is computationally very expensive. Fig. ^ shows 
Afcs and Afn for the same systems. Afes is calculated by explicitly evaluating the right hand 
side of (p^), whereas A/'n is obtained from extrapolating the transverse response from Fig. 
^ to g = 0. The difference between Af^s and Afn is Afs- For weak disorder {A/t < 5), the 
fluctuations from realization to realization are small. With increasing disorder, however, 
these fluctuations become quite large, as seen by the error bars in Fig. ^ which represent 
sample-to-sample fluctuations. Note that the zero-frequency transverse response shows little 
dependence on momentum in this approximation. 

An alternative definition of superfiuid density is as a stiffness to variations in the phase 
of the condensate wavefunction.0 Such a phase variation imposes a superfiuid velocity on 
the system, 

Vs = -Ve, (27) 

m 

where 6 is the phase of the condensate wavefunction. The total energy of the system increases 
due to the kinetic energy of the superflow, and is directly proportional to the density of 
superfiuid. The superfiuid density is then defined by 
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It is easy to estimate the superfluid density in this manner using the energy of the 
Hartree state (|^). (In principle, one should also include the change in the zero-point energy 
of the quasiparticles {\i7\), but this calculation is also costly.) To impose a twist in boundary 



conditions, we should change the hopping matrix elements tij to tije^^"-^ , where Aij is a 
vector potential whose sum along a path spanning the sample is 6. We should then solve the 
corresponding new non-linear Schrodinger equation, and compare the resulting condensate 
energies. Unfortunately, for 6 ^ mr this requires solving a complex non-linear Schrodinger 
equation. 

For 6 = 71, the Schrodinger equation remains real, but a new problem arises. Consider 
first the uniform case with V{i) = 0. With a vr phase twist, the ground state manifold of 
is doubly degenerate, and is spanned by the uniformly left- and right-moving condensates. 
This degeneracy frustrates our iterative convergence scheme, since linear combinations of 
these two degenerate solutions have spatially varying densities, driving even the Broyden 
method away from convergence. This problem persists in the disordered case. 

To avoid these complications, we note that for a phase difference of ^ = 27r the 
Schrodinger equation is unchanged. In the course of increasing the phase difference from 
to 27r, the ground state is deformed into the first excited state. Thus we have taken the 
energy ex=i of the first Hartree excited state in the absence of a twist to be the energy for 
introducing a 2tt phase twist across the sample. In our problem, the superfluid density can 
then be computed by 

2pst\^-^f. (29) 



(Strictly speaking, A/'es should be used in place of N in (p9D, since without disorder the nor- 
mal fluid density vanishes, while the longitudinal response is given by A/'cg. This correction 
is comparable in magnitude to the alteration of the zero-point motion of the quasiparticle 
energy, which we have also neglected in obtaining (|29|).) 

Fig. I shows the superfluid fraction obtained using ( PP] ) and again extrapolating to 
the thermodynamic limit by averaging over 7 realizations for system sizes from L = 8 to 
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18. For weak disorder, when the direct evaluation of the response function permits reliable 
extrapolation to q = 0, the two calculations agree. As the disorder grows, however, the 



A/t ~ 10 the direct calculation of Xij is no longer feasible because of the large number of 
samples required to obtain a reasonable statistical average. 

As shown in Fig. |^, systems with larger interaction are more robust to the addition of 
disorder and thus require larger A to reduce the superfluid density. (This was also the case 
with the condensate fraction; compare with Fig. ^.) In fact, the superfluid fraction remains 
substantial even when the Bogoliubov approximation begins to break down, i.e., when the 
condensate fraction becomes small. Note that the superfluid response involves both the 
bosons in the condensate and those that have been scattered out of it by interactions - 
when the condensate is accelerated, some of these scattered bosons accompany it. In the 
absence of disorder (and neglecting lattice effect^l) all bosons participate, and ps = p, 
even though the condensate can be substantially depleted. As Un/t is increased for fixed 
disorder, the condensate fraction A/q/A/" is reduced (because of increased scattering out of the 
condensate), while the superfluid fraction Ps/pis increased (because of decreased sensitivity 
to disorder when interactions are strong). 



At long wavelengths, the collective excitations of a uniform Bose condensate are phonons 
with a linear dispersion, uj = ck. For wavelengths comparable to the interparticle spacing, 
strongly interacting Bose fluids exhibit a roton minimum, and for even shorter wavelengths 
the collective excitations become ill defined, merging with the multiparticle continuum.0 

In the Bogoliubov approximation for lattice bosons without disorder, the quasiparticle 
spectrum can be solved analytically. The tight-binding dispersion is 



fluctuations in the response 




beyond 



VII. EXCITATION SPECTRUM 




(30) 



and quasiparticle dispersion is 
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= V2[/nek + . (31) 

A linear phonon dispersion holds for wavelengths that are long enough that (a) ek is much 
less than t, so that the tight-binding dispersion is nearly quadratic (ek ~ tk"^), and (b) ek is 
much less than 2Un, so that the first term in the square-root in (|3l|) dominates. The speed 



of sound c is then y/2Unt. 

The Bogoliubov approximation is too crude to capture the roton minimum found in 



real strongly-interacting Bose fluids, and for higher momenta the excitation energy (|3lD 
rises monotonically. At short wavelengths pair scattering can be neglected, and the quasi- 
particles behave as free particles with a Hartree energy uhf = tk"^ + Un. 

When disorder is introduced, translational invariance is destroyed, and momentum is no 
longer a good quantum numberJll Are the excitations created by 7^ localized? Even if the 
Hartree states (px are locahzed, the quasiparticle states created by 7^^ need not be - the 
condensate is extended, and can mediate non-local scattering via the inner product Sw in 
the Bogoliubov Hamiltonian ([T^). 

Since the quasi-particle operators 7^^ do not simply add a boson, but superpose a particle 
and a "hole" (a particle supplied by the condensate), the participation ratio used for the 
Hartree excitations (eq. (|)) is inappropriate. Transforming ([131) to the site basis, the 
quasiparticles are created by 



T.iU,^bl + V,A)■ (32) 



Adding an excitation in state fi to the ground state, there is amplitude f/^j to create a 
particle at site i and to create a "hole" there. The net particle density at site i in the 
one quasiparticle state 7^|G) differs from the density of the ground state itself by 

Sn^, = {G\j,nal\G) - {G\n,\G) 

= K^\' + \V^i\'. (33) 

The corresponding "participation ratio" specifying the degree of delocalization of this density 
fluctuation is then 
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\2 



m = wf^- (34) 



For an extended excitation, the participation ratio (^J) should scale linearly with the 
volume of the system; for a localized excitation, the participation ratio should become 
independent of the volume for systems larger than the localization length. Unfortunately, 
we could not perform a reliable scaling analysis with the small systems available to us, and 
we therefore could not infer the nature of the excitations created by 7^^. On general grounds, 
however, we expect the nature of the excitations in a disordered Bose condensate to be given 
by the localization problem for phononsS Thus in two dimensions, all excitations should 
be localized with a frequency dependent localization length .^(ti;) ~ exp[y4/ti;^] for arbitrary 
disorder. (In three dimensions, a mobility edge separates the extended low energy phonons 
from higher energy localized modes.) 

Within the Bogoliubov approximation (|10|), the excited states of the system are in- 
dependent bosons created by 7^^. If we assume that the temperature is low enough that 
thermal fluctuations do not change the excitation spectrum appreciably, but merely excite 
the quasiparticle states according to the Bose-Einstein distribution, the specific heat is then 
completely determined by the density of quasiparticle states per unit energy. For a density 
of states g{uj), the specific heat is 

C{T) g{u)e^l^BT^^ 

(35) 



kB Jo {keTy [e^l^BT _ i]2 • 

To infer the behavior of the specific heat at low temperatures, it is useful to introduce 
the "integrated density of states" 

N{uj) = r g{uj')du' = Y,Q{uj- u^), (36) 
■^0 fj, 

where Q{uj) is the Heaviside step function. N{uj) gives the number of states with energy less 
than or equal to u. As a monotonic function of uj, N{uj) is easier to fit than the spiky g{uj) 
for finite systems. If N{uj) ~ u^, then C{T) ~ for low temperatures. For a linear phonon 
dispersion u = ck in a (i-dimensional box of linear dimension L, the integrated density of 



states varies as N^u) ~ {Loj/cY, so that the specific heat of a uniform Bose condensate 
varies as T'^ at low temperatures. 

For the finite x Ly lattices we consider, the momenta k are restricted to a discrete set of 
allowed values. Only a limited number of these values satisfy the condition that ek be much 



smaller than both t and 2Un (or equivalently, that uj be much smaller than both y2Unt 
and 2Un) needed for ( pOD and (|31|) to yield the correct linear dispersion in the absence of 
disorder. To obtain enough states to permit a fit to the density of states in this regime, we 
are forced to work with large Un even though the Bogoliubov approximation is uncontrolled 
in this limit. The condition that ek is much less than t guarantees that we avoid the van 
Hove singularity at the center of the tight-binding band. (The van Hove singularity can 
also be pushed to higher energy by the judicious addition of further range hopping matrix 
elements which cancel the k'^ terms in (^) and prolong the k"^ dependence of Ck.) 

In the presence of disorder, the low energy integrated density of states will deviate from 
its pure uj'^ form. The integrated density of states divided by u"^ is shown for increasing 
disorder in Fig. ^ Each panel shows N{uj)/uP for one realization of disorder. The low- 
energy end of the spectrum (well below the van Hove singularity, and in the range which 
had a linear dispersion in the absence of disorder) is well-fit by 

N{u) =Au + 5cu^ (37) 
where A and B depend on both A/t and Un/t. 



Fig. |Ty shows the parameters A and B vs. disorder for Un/t = 3.3. (Qualitatively 
similar behavior is found for Un/t = 5.0 and 7.0.) For weak disorder the integrated density 
of states remains nearly quadratic in tu, consistent with the low energy excitations being 
weakly perturbed phonons. As disorder increases, however, a linear contribution to N{ijj) 
emerges, corresponding to a constant density of states g{uj). By the time the condensate is 
nearly completely depleted, the linear contribution to A^(ti;) dominates. 
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VIII. CONCLUDING REMARKS 



Gillis et alx^ have measured the low temperature (50 mK - IK) heat capacity of thin ^He 
films adsorbed in porous Vycor glass. At low coverages, they find that the heat capacity 
is linear, with no evidence of a superfiuid transition. This phase is identified with the 
insulating state of bosons localized by disorder, the "Bose glass." Above a critical coverage 
(corresponding to several monolayers), the low temperature phase is a superfiuid, with a 
heat capacity that varies as T^. 

Although the ^He is adsorbed as a few-monolayer film in the Vycor, the pores are 
connected to form a three-dimensional network. For sufficiently long wavelengths, three- 
dimensional behavior is expected. Why does the specific heat of the superfiuid vary as 

rather than T^, as expected for a three-dimensional condensate? The explanation is 
that for excitations with wavelengths less than the typical pore size Apore (several hundred 
Angstroms), the connectivity of the porous network is unimportant, and the density of 
states for phonon-like excitations will be that of a two-dimensional superfiuid. Thus above 
a crossover temperature ksT^ ~ he/ \ore-, the specific heat should vary as T^, until the 
roton contribution becomes appreciable. For an upper bound on T.j. we can use the bulk 
speed of sound c ~ 3 x 10'^ cm/sec, which gives a crossover temperature of 30 mK. This is 
surely an overestimate, since at the coverages studied by Gillis et al. the pores are not close 
to being filled and the compressibility is therefore much less than in bulk. An alternative 
estimate using the speed of sound in thin *^He films adsorbed on graphite gives a crossover 
temperature of 1 mK. 

How big a linear specific heat does one expect in the Bose glass? If one assumes a 
constant density of bosonic excitations in the Bose glass (as we found for the strongly 
disordered superfiuid), then the observed linear specific heat translates to roughly one mode 
per particle per 10 ^eV. This energy scale is comparable to the quantum confinement energy 
of a ^He atom trapped in a pore several hundred Angstroms in diameter. 

We have presented numerical solutions of the disordered Bose Hubbard model in the Bo- 
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goliubov approximation. This approximation correctly captures the long wavelength prop- 
erties of the clean Bose condensate, and is equivalent to the random phase approximation. 
It represents an expansion about the weakly disordered and weakly interacting limit, with 
a small parameter given by the depletion of the condensate. We find that weak disorder 
hardly affects the condensate fraction or the superfluid density. Instead, the condensate 
distorts to screen the imposed random potential. Interactions help stabilize the condensate, 
and prevent its collapse into a macroscopically occupied localized state. 

For strong disorder the condensate fraction and superfluid density are reduced, and 
ultimately vanish for sufficiently large disorder (although the Bogoliubov approximation is 
no longer controlled by this point). Our calculation therefore cannot access the critical 
properties of the superfluid-insulator transition. The Bogoliubov calculation, however, does 
suggest a promising variational state for the Bose glass. 

The clean Bose condensate has a linear low energy density of states in two dimensions, 
which implies a low temperature specific heat that varies as T^, as observed. We find that 
with increasing disorder a constant density of states appears at low energy. This constant 
density of states dominates as the condensate fraction and superfluid density become small, 
and leads to a linear low temperature specific heat. With our small sample sizes, we could 
not determine the extent to which these excitations are localized. 
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FIGURES 

FIG. 1. The condensate wavefunction 0o is concentrated at the minima of the "bare" poten- 
tial V{i). The effective potential felt by the bosons is therefore increased in the regions of high 
condensate density, so that the resulting self-consistently determined potential V{i) resembles the 
original potential V{i) with its lowest values lopped off. (Here Un/t = 3.3 and A/t = 10.0) 

FIG. 2. The autocorrelation function of the screened potential. Note that the screened 
potential does not develop any long-range correlations. (V = 210 sites and Un/t = 3.3) 

FIG. 3. An example of the participation ratio for the Hartree quasiparticle states (j)x as a 
function of their energy for several system sizes. {Un/t = 3.3 and A/t = 10.0) 

FIG. 4. The condensate fraction Nq/N versus disorder A/t for several interaction strengths, 
extrapolated to V — ^ oo. 

FIG. 5. The density {G\ni\G) of the Bogoliubov ground state versus the density A/'|^o(OP 
of the Hartree state. Note that despite the large depletion of the condensate in the Bogoliubov 
state, the total density is well approximated by the density of the completely condensed state. 
{Un/t = 3.3) 

FIG. 6. The transverse component of the current-current response function for the 12 x 13 
lattice. {Un/t = 3.3) 

FIG. 7. Mes/M and Mn/M vs. disorder for the 12 x 13 lattice. Ms = M^s - Mn- {Un/t = 3.3) 

FIG. 8. The superfiuid fraction as obtained through the twist method for several interaction 
strengths, extrapolated to V — oo. 

FIG. 9. N{uj)/uj^ vs. u for Un/t = 3.3 and {a)A/t = 0.0 {h)A/t = 10.0 {c)A/t = 18.0. As 
disorder increases, N{ijo)/uj^ diverges, indicating a deviation from the form N{u}) ~ u)'^. (V = 210 
sites.) 
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FIG. 10. Coefficients of the linear and quadratic parts of the integrated density of states 
N{u}) = Alo + BuP' vs. disorder (210 site system). As disorder is increased, the Unear term 
develops, indicative of a glassy system. 
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